/* 
Table 1: Estimates of the Effect of Veteran Status on Homeownership by Race (1960)
Figure 4: IV Estimates of Effect of Home Loan Guaranty on Probability of Homeownership (1960)

Chinyere O Agbai 

*/
matrix drop _all
macro drop _all

global dir = "YourFilePath"
global output = "$dir"
global working_data = "$dir/ipums60_pub.dta"

cd "$output"
use "$working_data", clear

set scheme plotplain

* State Fixed Effects Global 
global state_fe state_cat2 state_cat3 state_cat4 state_cat5 state_cat6 ///
	state_cat7 state_cat8 state_cat9 state_cat10 state_cat11 state_cat12 state_cat13 ///
	state_cat14 state_cat15 state_cat16 state_cat17 state_cat18 state_cat19 state_cat20 ///
	state_cat21 state_cat22 state_cat23 state_cat24 state_cat25 state_cat26 state_cat27 ///
	state_cat28 state_cat29 state_cat30 state_cat31 state_cat32 state_cat33 state_cat34 /// 
	state_cat35 state_cat36 state_cat37 state_cat38 state_cat39 state_cat40 state_cat41 ///
	state_cat42 state_cat43 state_cat44 state_cat45 state_cat46 state_cat47 state_cat48 ///
	state_cat49 state_cat50

***************************************************************************************
* Table 1: Estimates of the Effect of Veteran Status on Homeownership by Race (1960) * 
***************************************************************************************
* table headers 
	set more off
	putexcel set T1.xls, replace
	putexcel B2 = ("Outcome(Y)") C2 = ("All")  D2 = ("White")  E2 = ("Black")	
	putexcel A3 = ("Stage 1") A5 = ("Reduced Form") A7 = ("IV") A9 = ("OLS") A11 = ("N") ///
		B3 = ("Veteran Status") B5 = ("Homeownership") B7 = ("Homeownership") B9 = ("Homeownership")
	putexcel A1 = ("Estimates of the Effect of Veteran Status on Homeownership by Race, 1912-1932 Birth Cohorts")	
	putexcel A12 = ("Robust standard errors in parentheses. Models include state fixed effects.") 
	putexcel A13 = ("Source: 1960 IPUMS 1% Sample")
	putexcel A14 = ("*** p<0.001, ** p<0.01, * p<0.05, + p<0.10")
	
	
* STAGE 1  *********************************************************************
	* Full sample * 
	reg vetbin_60 inst_1927 age $state_fe [weight = perwt], vce(cluster statefip)	
		matrix list r(table)
		matrix define reg_output= r(table) 
		local b = string(reg_output[1,1], "%04.3f") 
		local pval =  string(reg_output[4,1], "%04.3f") 
		
		local pval_star = ""
			if (`pval' < 0.001) local pval_star = "***"
			else if (`pval' < 0.01) local pval_star = "**"
			else if (`pval' < 0.05) local pval_star = "*"
			else if (`pval' < 0.10) local pval_star = "+"

		local row = 3
		putexcel C`row'= ("`b'`pval_star'")   
		
		local row = `row' + 1  
			
		local se = string(reg_output[2,1], "%04.3f")
		putexcel C`row'= ("(`se')") 
				
	* WHITE * 
	reg vetbin_60  inst_1927 age $state_fe [weight = perwt] if race ==1, vce(cluster statefip)	
	di "race == 1"
	
		matrix list r(table)
		matrix define reg_output= r(table) 
		local b = string(reg_output[1,1], "%04.3f") 
		local pval =  string(reg_output[4,1], "%04.3f") 
		
		local pval_star = ""
			if (`pval' < 0.001) local pval_star = "***"
			else if (`pval' < 0.01) local pval_star = "**"
			else if (`pval' < 0.05) local pval_star = "*"
			else if (`pval' < 0.10) local pval_star = "+"
		
		local row = 3
		putexcel D`row'= ("`b'`pval_star'")   
		
		local row = `row' + 1  
			
		
		local se = string(reg_output[2,1], "%04.3f")
		di `se'
		putexcel D`row'= ("(`se')") 
		
	* Black * 
	reg vetbin_60  inst_1927 age $state_fe [weight = perwt] if race == 2, vce(cluster statefip)	
	
		matrix list r(table)
		matrix define reg_output= r(table) 
		local b = string(reg_output[1,1], "%04.3f") 
		local pval =  string(reg_output[4,1], "%04.3f") 
		
		local pval_star = ""
			if (`pval' < 0.001) local pval_star = "***"
			else if (`pval' < 0.01) local pval_star = "**"
			else if (`pval' < 0.05) local pval_star = "*"
			else if (`pval' < 0.10) local pval_star = "+"
		
		local row = 3
		putexcel E`row'= ("`b'`pval_star'")   
		
		local row = `row' + 1  
		
		local se = string(reg_output[2,1], "%04.3f")
		putexcel E`row'= ("(`se')") 
		

* REDUCED FORM 1  *********************************************************************
	* FULL SAMPLE * 
	reg ownbin  inst_1927 age $state_fe [weight = perwt], vce(cluster statefip)
		matrix list r(table)
		matrix define reg_output= r(table) 
		local b = string(reg_output[1,1], "%04.3f") 
		
		di `b'
		local pval =  string(reg_output[4,1], "%04.3f") 
		di `pval'
		
	local pval_star = ""
		if (`pval' < 0.001) local pval_star = "***"
		else if (`pval' < 0.01) local pval_star = "**"
		else if (`pval' < 0.05) local pval_star = "*"
		else if (`pval' < 0.10) local pval_star = "+"
		
	local row = 5
		
	putexcel C`row'= ("`b'`pval_star'")   
		
	local row = `row' + 1  	
		
	local se = string(reg_output[2,1], "%04.3f")
	putexcel C`row'= ("(`se')") 
	
	* WHITE * 
	reg ownbin  inst_1927 age $state_fe [weight = perwt] if race == 1, vce(cluster statefip)
	di "race == 1"
	
		matrix list r(table)
		matrix define reg_output= r(table)
		local b = string(reg_output[1,1], "%04.3f") 
		di `b'
		local pval =  string(reg_output[4,1], "%04.3f") 
		di `pval'
		
	local pval_star = ""
		if (`pval' < 0.001) local pval_star = "***"
		else if (`pval' < 0.01) local pval_star = "**"
		else if (`pval' < 0.05) local pval_star = "*"
		else if (`pval' < 0.10) local pval_star = "+"
	
	local row = 5
		
	putexcel D`row'= ("`b'`pval_star'")   
		
	local row = `row' + 1  	
		
	local se = string(reg_output[2,1], "%04.3f")
	di `se'
	putexcel D`row'= ("(`se')") 

	* BLACK * 
	reg ownbin  inst_1927 age $state_fe [weight = perwt] if race == 2, vce(cluster statefip)
	di "race == 2"
	
		matrix list r(table)
		matrix define reg_output= r(table) 
		local b = string(reg_output[1,1], "%04.3f") 
		di `b'
		local pval =  string(reg_output[4,1], "%04.3f") 
		di `pval'
		
	local pval_star = ""
		if (`pval' < 0.001) local pval_star = "***"
		else if (`pval' < 0.01) local pval_star = "**"
		else if (`pval' < 0.05) local pval_star = "*"
		else if (`pval' < 0.10) local pval_star = "+"
	
	local row = 5
		
	putexcel E`row'= ("`b'`pval_star'")   
		
	local row = `row' + 1  	
		
	local se = string(reg_output[2,1], "%04.3f")
	di `se'
	putexcel E`row'= ("(`se')") 


* S2 (Fuzzy Regression Discontinuity)  *****************************************
	* FULL SAMPLE * 
	ivreg ownbin age $state_fe (vetbin_60= inst_1927) [weight = perwt], cluster (statefip)
	predict s2_all_60 if e(sample)
	predict se_s2_all_60 if e(sample), stdp 
				
	local pval = 2 * ttail(e(df_r), abs(_b[vetbin] / _se[vetbin]))
	local pval_star = ""
		if (`pval' < 0.001) local pval_star = "***"
		else if (`pval' < 0.01) local pval_star = "**"
		else if (`pval' < 0.05) local pval_star = "*"
		else if (`pval' < 0.10) local pval_star = "+"
		
	local row = 7 
		
	local b = string(_b[vetbin_`i'] ,"%12.3f")
	putexcel C`row'= ("`b'`pval_star'")   
			
	local row = `row'+1  
			
	local se = string(_se[vetbin_`i'], "%12.3f")
	putexcel C`row'= ("(`se')") 
	

	* WHITE * 
	ivreg ownbin age $state_fe (vetbin_60= inst_1927) [weight = perwt] if race == 1, cluster (statefip)
	predict s2_white_60 if e(sample)
	predict se_s2_white_60 if e(sample), stdp 
	
	local pval = 2 * ttail(e(df_r), abs(_b[vetbin] / _se[vetbin]))
	local pval_star = ""
		if (`pval' < 0.001) local pval_star = "***"
		else if (`pval' < 0.01) local pval_star = "**"
		else if (`pval' < 0.05) local pval_star = "*"
		else if (`pval' < 0.10) local pval_star = "+"
			
	local row = 7 

	local b = string(_b[vetbin_`i'] ,"%12.3f")
	putexcel D`row'= ("`b'`pval_star'")   
			
	local row = `row'+1  
			
	local se = string(_se[vetbin_`i'], "%12.3f")
	putexcel D`row'= ("(`se')") 
		
		
	
	* BLACK * 
	ivreg ownbin age $state_fe (vetbin_60= inst_1927) [weight = perwt] if race == 2, cluster (statefip)
	predict s2_black_60 if e(sample)
	predict se_s2_black_60 if e(sample), stdp 
	
	local pval = 2 * ttail(e(df_r), abs(_b[vetbin] / _se[vetbin]))
	local pval_star = ""
		if (`pval' <= 0.001) local pval_star = "***"
		else if (`pval' <= 0.01) local pval_star = "**"
		else if (`pval' <= 0.05) local pval_star = "*"
		else if (`pval' <= 0.10) local pval_star = "+"
		
	local row = 7 
		
	local b = string(_b[vetbin_`i'] ,"%12.3f")
	putexcel E`row'= ("`b'`pval_star'")   
			
	local row = `row'+1  
			
	local se = string(_se[vetbin_`i'], "%12.3f")
	putexcel E`row'= ("(`se')") 

* OLS **************************************************************************
	* FULL SAMPLE * 	
	reg ownbin vetbin_60 age $state_fe [weight = perwt] , vce(cluster statefip)
		matrix list r(table)
		matrix define reg_output= r(table)
		local b = string(reg_output[1,1], "%04.3f") 
		local pval =  string(reg_output[4,1], "%04.3f") 
		local obs = string(e(N), "%12.0gc")
		di "b = `b'; pval = `pval'; N=`obs'"


	local pval_star = ""
		if (`pval' < 0.001) local pval_star = "***"
		else if (`pval' < 0.01) local pval_star = "**"
		else if (`pval' < 0.05) local pval_star = "*"
		else if (`pval' < 0.10) local pval_star = "+"
		
	local row = 9
		
	putexcel C`row'= ("`b'`pval_star'")   
		
	local row = `row' + 1  	
		
	local se = string(reg_output[2,1], "%04.3f")
	putexcel C`row'= ("(`se')") 
	
	local row = `row' + 1  	
	putexcel C`row'= ("`obs'") 

	
	* WHITE * 	
	reg ownbin vetbin_60 age $state_fe [weight = perwt] if race == 1, vce(cluster statefip)
	
		matrix list r(table)
		matrix define reg_output= r(table) 
		local b = string(reg_output[1,1], "%04.3f") 
		local pval =  string(reg_output[4,1], "%04.3f") 
		local obs = string(e(N), "%12.0gc")
		di "b = `b'; pval = `pval'; N=`obs'"
		
	local pval_star = ""
		if (`pval' < 0.001) local pval_star = "***"
		else if (`pval' < 0.01) local pval_star = "**"
		else if (`pval' < 0.05) local pval_star = "*"
		else if (`pval' < 0.10) local pval_star = "+"
	
	local row = 9
		
	putexcel D`row'= ("`b'`pval_star'")   
		
	local row = `row' + 1  	
		
	local se = string(reg_output[2,1], "%04.3f")
	di `se'
	putexcel D`row'= ("(`se')") 
	
	local row = `row' + 1  	
	putexcel D`row'= ("`obs'") 

	* BLACK * 	
	reg ownbin vetbin_60 age $state_fe [weight = perwt] if race == 2, vce(cluster statefip)
	
		matrix list r(table)
		matrix define reg_output= r(table) 

		local b = string(reg_output[1,1], "%04.3f") 
		local pval =  string(reg_output[4,1], "%04.3f") 
		local obs = string(e(N), "%12.0gc")
		di "b = `b'; pval = `pval'; N=`obs'"
		
	local pval_star = ""
		if (`pval' < 0.001) local pval_star = "***"
		else if (`pval' < 0.01) local pval_star = "**"
		else if (`pval' < 0.05) local pval_star = "*"
		else if (`pval' < 0.10) local pval_star = "+"
	
	local row = 9
		
	putexcel E`row'= ("`b'`pval_star'")   
		
	local row = `row' + 1  	
		
	local se = string(reg_output[2,1], "%04.3f")
	putexcel E`row'= ("(`se')") 
	local row = `row' + 1  	
	putexcel E`row'= ("`obs'") 


****************************************************************************************************
* Figure 4: : IV Estimates of Effect of Home Loan Guaranty on Probability of Homeownership (1960) * 
****************************************************************************************************

* xlabel for race and veteran status 1960 
	gen racevet = .
	replace racevet = 1 if vetbin_60 == 1 & race == 1
	replace racevet = 2 if vetbin_60 == 1 & race == 2
	replace racevet = 3 if vetbin_60 == 0 & race == 1
	replace racevet = 4 if vetbin_60 == 0 & race == 2
			
	label define racevet_cat1 1 "White Vets 1960" 2 "Black Vets 1960" 3 "White Non-Vets 1960" 4 "Black Non-Vets 1960"
	label values racevet racevet_cat1
	label var racevet "Race & Veteran Status"
		
gen racevet_mean = .
gen lb_pi = .
gen ub_pi = .

forvalues i = 1/4 {
	sum s2_white_60 if racevet == `i' 
	replace racevet_mean = r(mean) if racevet == `i' & racevet_mean ==.
	sum s2_black_60 if racevet == `i'
	replace racevet_mean = r(mean) if racevet == `i' & racevet_mean ==.
	di "racevet cat = `i' & mean = `r(mean)'"
	
	replace lb_pi = s2_white_60 - invttail(e(df_r), 0.025) * se_s2_white_60 if racevet == `i' & lb_pi ==.
	replace ub_pi = s2_white_60 + invttail(e(df_r), 0.025) * se_s2_white_60 if racevet == `i' & ub_pi ==.

	replace lb_pi = s2_black_60 - invttail(e(df_r), 0.025) * se_s2_black_60 if racevet == `i' & lb_pi ==.
	replace ub_pi = s2_black_60 + invttail(e(df_r), 0.025) * se_s2_black_60 if racevet == `i' & ub_pi ==.
	
	sum lb_pi if racevet == `i'  
	replace lb_pi = r(mean) if racevet == `i' 
	sum ub_pi if racevet == `i'
	replace ub_pi  = r(mean) if racevet == `i' 
	}
	format racevet_mean %8.2f 

graph twoway (bar racevet_mean racevet if racevet==1, color(eltgreen) barw(0.9)) ///
	(rcap lb_pi ub_pi racevet if racevet ==1,  vertical lcolor(gs4)) ///
	(bar racevet_mean racevet if racevet==2, color(lavender) barw(0.9)) ///
	(rcap lb_pi ub_pi racevet if racevet ==2,  vertical lcolor(gs4)) ///
	(bar racevet_mean racevet if racevet==3, color(eltgreen) barw(0.9)) ///
	(rcap lb_pi ub_pi racevet if racevet ==3,  vertical lcolor(gs4)) ///
	(bar racevet_mean racevet if racevet==4, color(lavender) barw(0.9)) ///
	(rcap lb_pi ub_pi racevet if racevet ==4,  vertical lcolor(gs4)) ///
	(scatter racevet_mean racevet, mlabel(racevet_mean) msymbol(none) mlabposition(1)), ///
	ylabel(0(0.1)0.8, labsize(medium)) ytitle("Predicted Pr(Homeownership)", size(medium)) ///
	xlabel(1.5 "Veteran" 3.5 "Non-Veteran", noticks labsize(medium)) xtitle("") legend(off)

graph export "fig4.png", //replace

*  age & age2 models 
	ivreg2 ownbin age $state_fe (vetbin_60= inst_1927) [weight = perwt], cluster (statefip)
	estat ic //aic 313322
	weakivtest 
	
	ivreg2 ownbin age age2 $state_fe (vetbin_60= inst_1927) [weight = perwt], cluster (statefip)
	estat ic //aic 315487

* Weak IV test
ivregress 2sls ownbin age $state_fe (vetbin_60= inst_1927) [weight = perwt], cluster (statefip)
estat firststage
